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We describe in detail our calculation of the full supersymmetric (SUSY) QCD corrections to 
neutralino annihilation into heavy quarks and extend our numerical analysis of the resulting dark 
matter relic density to scenarios without scalar or gaugino mass unification. In these scenarios, the 
W), final state is often composed of top quarks and the annihilation proceeds through Z°-boson or scalar 

top-quark exchanges. The impact of the corrections is again shown to be sizable, so that they must 
' be taken into account systematically in global analyses of the supersymmetry parameter space. 
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I. INTRODUCTION 



The search for physics beyond the Standard Model (SM) is no longer restricted to colliders only. In fact, the most 
compelling evidence for new physics comes today from cosmological observations such as the mission of the Wilkinson 
Microwave Anisotropy Probe (WMAP), which have determined the matter and energy decomposition of our universe 
with unprecedented precision. These observations indicate the existence of Cold Dark Matter (CDM) in the universe, 
\ which cannot be accounted for by the SM and likely consists of Weakly Interacting Massive Particles (WIMPs) with 
non-relativistic velocities. A combination of the five-year measurement of the cosmic microwave background by the 
WMAP mission with supernova and baryonic acoustic oscillation data yields the narrow 2c-interval for the relic 
C*~) density of dark matter [Ij 
O ■ 

O; 0.1097 < n CDM h 2 < 0.1165, (1) 

where h denotes the present Hubble expansion rate Hq in units of 100 km s 1 Mpc 1 . 
Q\ \ The measured relic density of dark matter can be used to constrain extensions of the Standard Model, which provide 
a viable WIMP candidate. In the Minimal Supersymmetric Standard Model (MSSM) with imparity conservation, this 
could be the Lightest Supersymmetric Particle (LSP) x, if it is neutral and a color singlet. One can then calculate 
its relic density, compare it with the experimental limits in Eq. ([1]), and identify the favored regions of the MSSM 
parameter space. The relic density 

n , h 2 = W*. (2) 

Pc 

is proportional to the present number density uq and the mass of the LSP. p c = 3Hq / (8irGtq) is the critical 
density of our universe, and Gn is the gravitational constant. The present number density no is obtained by solving 
the Boltzmann equation describing the time evolution of the number density 



driy 
~dt 



v = - 3Hn^ - (a ann v) (n\ - n^ q ) . (3) 



The first term on the right-hand side corresponds to a dilution due to the expansion of the universe, and the second 
term corresponds to a decrease due to annihilations and co-annihilations of the relic particle into SM particles 0. 
Here, H denotes the time-dependent Hubble expansion parameter, and n eq the density of the relic particle in thermal 
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FIG. 1: Tree-level Feynman diagrams for the annihilation of a neutralino pair into a quark-antiquark pair through the exchange 
of a Z°-boson, a neutral Higgs-boson = (h°, H°, A ), or a squark q~i (i — 1,2). 

equilibrium. Details of the dark matter interactions enter the Boltzmann equation through the thermally averaged 
cross section (cr ann i>). The cross section takes into account the thermal velocity distribution of the relic particle and 
is calculated for a given temperature T by 



where K\ and K% are the modified Bessel functions of the first and second kind, respectively. The center-of-momentum 
energy y/s is related to the particle mass and the relative momentum p cm of the annihilating pair through 



In order to keep up with current and future experimental improvements, one has to understand and reduce the 
different uncertainties involved in the analysis, both for the prediction of the dark matter relic density and the 
extraction of new mass parameters from cosmological data. These uncertainties include, e.g., a modification of the 
Hubble expansion rate due to quintessence or an effective dark energy density [H, differences in the new physical 
particle masses obtained with different spectrum codes Q, or a lack of precision in the annihilation cross-section 
of dark matter particles [f|. In this paper, we focus on the impact of next-to- leading order QCD and SUSY-QCD 
corrections on the latter, but other possible uncertainties will also be briefly discussed. 

In many scenarios of the MSSM, the lightest neutralino is the LSP and therefore a suitable dark matter candidate. 
The thermally averaged cross section is then obtained by computing all relevant neutralino annihilation cross sections 
into SM particles. Most prominent arc the processes with two-particle final states such as a fermion-antifermion 
pair or a combination of gauge (W ± ,Z°) and Higgs bosons (h°, H°, A°,H ± ) [1,0. In this paper, we focus on the 
annihilation into a massive quark-antiquark pair, since the leading order cross section with a fermion-antifermion final 
state is proportional to the mass of the fermion, which disfavors the light quarks. Moreover, the annihilation into 
heavy quarks is important in the regions of parameter space allowed by Eq. ([!} . We now present the full details of 
our calculation and investigate scenarios with dominant top-quark final state contributions, extending our analysis 
of Refs. [1, [^] beyond minimal supergravity (mSUGRA) models. In mSUGRA models, one is constrained by having 
only five universal high-scale parameters [mo, rn 1 / 2 , ^4oj tan/3, and sgn(/x)], and in regions of parameter space, where 
quark final states are important, the cross section is dominated by an exchange of Higgs bosons. Here, we relax 
the unification of either the scalar masses or the gaugino masses, one at a time. This allows for scenarios different 
from mSUGRA, where the annihilation cross section is not necessarily dominated by Higgs-boson exchanges. Apart 
from a Higgs-boson dominated scenario, we thus analyze scenarios, where Z°-boson exchanges in the s-channel or 
squark exchanges in the t- and u-channels play an important role. The full QCD and SUSY-QCD corrections in these 
scenarios turn out to be significant, and we have therefore included them into the public code micrOMEGAs [lfj. 

This paper is organized as follows: In Sec. Inland appendices [Al and [Bl we give all necessary details of the calculation 
related to the virtual loop corrections, the renormalization procedure, and the calculation of real gluon emission, in 
particular the subtraction of the induced soft and collinear singularities. We then continue in Sec. IIIII with a discussion 
of the MSSM models beyond scalar or gaugino mass unification. In Sec. IIV1 we analyze the impact of the radiative 
corrections on the relic density in these models. Finally, our results are summarized inlVl 



The annihilation of neutralinos into quarks proceeds at tree-level through an exchange of a Z -boson and Higgs- 
bosons in the s-channel and through the exchanges of scalar quarks in the t- and u-channels (see Fig.QJ. By definition, 
a WIMP can only have electroweak interactions, and in order to reach a sufficient annihilation rate to explain the 
dark matter relic density, the cross section has to be enhanced, e.g., by a resonance. One often finds that there is 




(4) 



S = 4(m|+^ ro ) Hi- 
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FIG. 2: SUSY-QCD loop diagrams contributing to the annihilation of neutralinos into quarks. 



one contribution which dominates the whole cross section. This allows us, by choosing different scenarios, to study 
contributions from various channels. Moreover, it allows us also to isolate effects of radiative corrections coming from 
different sources. 

We have computed the full QCD and SUSY-QCD corrections to neutralino annihilation into quarks. The next-to- 
lcading order cross section contains virtual contributions stemming from loop diagrams and real contributions, which 
are due to the radiation of an additional gluon. Symbolically, the next-to-leading order (NLO) cross section can be 
written as 

a NLO =f da v +f da R , (5) 
where cr V ' R denote the virtual and real emission parts integrated over the two- and three-particle phase space, respec- 
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FIG. 3: SUSY-QCD counterterm diagrams contributing to the annihilation of neutralinos into quarks. 



tively. The cross section er v , which includes the tree- level and the one- loop virtual corrections, is given by 



.00 



1 



K{s,m 2 q ,m 2 q ) 
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(0) 



where k(x, y, z) — yjx — y — z) 2 — Ayz. For our notation and conventions, we refer the reader to App. [X] There are 
two types of contributions to the one- loop amplitude A^x-loopj those coming from the loop diagrams and those coming 
from the counterterms. The loop diagrams, depicted in Fig. [2j contain ultraviolet (UV) and infrared (IR) divergent 
loop integrals. We regulate both types of divergences dimensionally (d = 4 — 2e) and evaluate the loop integrals in 
the dimensional reduction scheme (DR). The full analytic results for the loop diagrams are given in App. [B] 

The UV divergences are compensated by counterterms (see Fig. [3]) related to quarks and their scalar super-partners, 
the squarks. All counterterm vertices with quark or squark legs contain wave-function renormalization factors SZ^, 
6Zr and SZij, that follow from replacing the (s)quark fields by 
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qL 
qn 



q, 



(S i:j + -5Zij)qj . 



(J) 



Although in principle only the wave-function renormalization constants of the external quarks have to be included, we 
also include the squark renormalization constants, as they allow us to perform simpler UV-convergence checks. In the 
full calculation, the squark wave-function renormalization constants cancel out. The wave-function renormalization 
constants are determined by requiring the residues of the propagators to remain at unity even at one-loop order. This 
condition gives 
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(8) 
(9) 

(10) 



where HL,n(k 2 ) and ^lsL,SR(k 2 ) stand for the vector and the scalar parts of the two-point Green's function as defined 
in Ref. [ll| and tl(m 2 ) — [^l^n^fc 2 )] fe2 _ m 2- After performing the wave-function renormalization, the remaining 
divergences are canceled by renormalizing the coupling constants. In our case, the coupling counterterms that receive 
contributions proportional to the strong coupling constant a s are the quark Yukawa couplings through the masses 
of the quarks, the squarks masses and the squark mixing angle 9 q . A very important contribution, in particular in 
scenarios with a dominant Higgs-boson exchange, comes from renormalizing the Yukawa couplings of the quarks 0]. 
In our calculation, we use the DR Yukawa couplings for both the top and the bottom quarks. As the quark masses 
that serve as inputs are defined in different schemes, we take two distinct approaches for top and bottom quarks. For 
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top quarks, the input is the on-shell mass mt = 172.4 GeV measured at the Tevatron |12j, and the DR-mass of the 
top quark (and hence the DR Yukawa coupling) is obtained by subtracting the finite on-shell countertcrm 



Sm = —'St 
9 2 



m q (U L (m q ) + Il R (m q )) + U S L(rn q ) + U S R(m q 



(11) 



On the other hand, the input mass mt,(mb) for bottom quarks is extracted in the MS renormalization scheme from 
the Standard Model analysis of T sum rules [l3j]. In order to obtain the appropriate bottom Yukawa coupling in 
the DR renormalization scheme within the MSSM, we first use the Standard Model next-to-next-to-leading order 
(NNLO) renormalization group evolution to obtain the mass of the bottom quark at the scale Q — 2m^ |14|. Still 

in the SM, we then convert m b ls (Q) to m b JB -(Q) [14], and finally we apply the threshold corrections including also 
contributions from SUSY particles in the loop. For the last step, we take into account the fact that the sbottom-gluino 
and stop-chargino one-loop contributions are considerably enhanced for large tan (3 or large Af, and can be resummed 
to all orders in perturbation theory [TEl IT(| . Denoting the resummable part by and the finite one-loop remainder 
by Am;,, the bottom quark mass is then given by 



DR, MSSM 



DR, SM 



(Q) 



(Q) 



1 + A b 



In the squark sector, we pick five independent quantities, m^, m^ 2 , mg 
symmetry. We renormalize the masses of three squarks in the on-shell scheme, which leads to the counterterm 



Am b . (12) 
r , 9z, and Or, in order to respect the SU(2) 



5mi. = St 



n|(m|) 



(13) 



The remaining mass of the heavier scalar bottom quark mg 2 is treated as dependent, 
1 



6m? cos 2 6 1 + Sm^ sin 2 % - 5m| cos 2 9- b + (m~ 2 - m| ) sin 2% S9t - (m^ - m\ ) sin 29 b 59 b . (14) 



The squark mixing angle is renormalized in the DR scheme and so the corresponding counterterms of the squark 
mixing matrices Rf, contain only the divergent parts. The counterterms can be determined as 



(15) 



fc=i 



using only the divergent parts of the wave-function renormalization constants. This is equivalent to fixing the mixing 
angle as 
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(16) 



After renormalization, the virtual cross section cr v still contains IR divergences, that come from an exchange of a 
gluon in the loops. In order to compensate them, one has to include the real cross section er R coming from diagrams 
with an additional gluon in the final state (see Fig. [4j . The symbolic formula of Eq. ([5]) cannot be directly applied, 
since the divergences appear as 1/e and 1/e 2 poles in the number of dimensions in the one-loop amplitude and as a 
divergent matrix element in certain regions of the 2 — ► 3 parameter space. A convenient way to combine these two 
cross sections and to cancel the divergences is the dipole subtraction method [T^. Using this method, we can write 
the total next-to-leading order cross section cr NLO as 
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(17) 



where we introduced an auxiliary cross section a A . The auxiliary cross section does not contribute to the total cross 
section and serves only as a tool to cancel the IR divergences. It has the same divergence structure as the real cross 
section and at the same time its structure allows for a partial integration of the gluon phase space, so that it can also 
be added to the virtual cross section, canceling the divergences in each part. The 2 — > 3 matrix element leading to 
the auxiliary cross section for the real part is constructed from two dipole contributions X>3i % and D32 i as 



\M 



2— 3|2 
aux I 



£>3i,2(fci, k 2 , k 3 ) + V 32t i{ki,k 2 , k 3 ) . 



(18) 



FIG. 4: Bremsstrahlung diagrams contributing to the annihilation of neutralinos into quarks. 



where the fc, are the four-momenta of the final state quarks and of the gluon. The dipole contributions 1531,2 and 
D32 i are related by a simple interchange k\ <-> k 2 with P31 2 given by 



^31,2(^1,^2, k 3 ) = C F 



8na s 
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2(1 - 2 M 



1 — X2 2 — XI — X2 



4 - 4 M 2 1 _ 2m 2 



2/y, 



iri — 1 
x 2 — 2\i 2 1 — x 2 



(19) 



where aij = 2ki.q/q 2 , (j, q = m q /*/s, and q 2 = (fci + k 2 + k 3 



appearing in Eq. (|19[) is calculated using different kinematics with redefined 4-momcnta 
1 



s [17|. The leading order matrix element \Ait 



K l ^ ""31 



1 - 4 M 2 



1 



2 - 
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1 - 4 M 2 
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/to 



2 - 2 X2 ^ 



The auxiliary matrix element that cancels the infrared divergences of the virtual matrix element is 



1*- 1 - 2C *£w±7)\ M * 



S12 



V q (si2,m q ,m q ;e) 



1 , . 3 jJj 2 it 2 

+ —r q ( mq ;e) + -\n— + 5-- 



where s%2 — s — 2m 2 . The function V q is composed of a singular part, 



V (s) (si2,m 9 ,TO,j;e) = 



l+(3 2 
2/3 



1,1-/3 1 , 2 1 - P tt 2 1-0 2 
— In In h in In tt 

e 1 + 2 1 + 6 1 + l + P 2 



with = ^Jl — Afi 2 , and a non-singular part given by 



V^ s H Sl2 ,m qi m q ) = -In 
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+ ln(l - -V^^) - 21n(l - Vl^P) - - 



1-/3 2 



In 



V^W 2 



1 + /3 2 2 



The function r 9 is defined as 



T q (m q ; e) = Cf 



f f 

7 + 2 ln #" 2 



(23) 



(24) 



where Q is the renormalization scale. Apart from convergence checks of our calculation, we also performed checks of 
the finite part of the one-loop diagrams against the results obtained by automatic computer packages FeynArts and 
FormCalc [l8j] and parts of the calculation also against existing results, e.g. in Ref. (T^. The UV and IR divergent 
parts of our loop integrals were checked and agree with the results given in Ref. (2fJ . 



III. SUSY MODELS BEYOND MINIMAL SUPERGRAVITY 



In our previous publications, we studied the impact of the SUSY-QCD corrections to neutralino annihilation within 
minimal supergravity (mSUGRA) models, defined at the grand unification scale A/gut in terms of a universal scalar 
mass mo, a universal gaugino mass m.1/2, a common trilinear coupling Aq, the ratio tan/3 of the Higgs doublet vacuum 
expectation values, and the sign of the higgsino mass parameter \i In mSUGRA models, the lightest neutralino 

is often a P-ino and annihilates preferably through resonant Higgs-boson exchanges. For example, in the focus-point 
region at high scalar masses Too the cross section is dominated by heavy CP-even Higgs bosons decaying into top 
quarks. Regions with important bottom quark final states include those with small gaugino masses TO4/2, where 
light CP-even Higgs-boson exchanges dominate, and the A°-funnel region at high values of tan (3, where the bottom 
Yukawa coupling is enhanced and the annihilation proceeds through pseudo-scalar Higgs-boson exchanges. At large 
tan (3, the bottom-quark contribution can also become sizable in the focus-point region. 

In this work, we study in detail numerically the annihilation of neutralinos into top quark-antiquark pairs through 
the exchange of Z°-bosons in the s-channel and of top squarks in the t- and it-channels (see Fig. [IJ. These channels 
can be enhanced in models, where some of the unification conditions have been relaxed, which is well motivated 
theoretically [2l], [2^, [23j] . To be concrete, we focus on two sets of models based on the SO(10) Grand Unified Theory 
(GUT): models with non-universal Higgs masses (NUHM), and models without gaugino mass unification. SO(IO) 
theories are particularly promising, as they involve complete 16-dimensional matter multi plet s with a right-handed 
neutrino and can be embedded in string theories involving larger groups like E 8 or SO(32) tHIH- 



A. SUSY models with non-universal Higgs masses 

In SO(IO) SUSY GUTs, the matter superfields of one generation belonging to a 16-dimensional representation are 
completely mass degenerate, if the SUSY-breaking masses are acquired above the SO(I0)-breaking scale. Flavor-blind 
mechanisms can furthermore lead to universal masses of all matter scalars. However, if the Higgs doublets H u and H d 
belong to different, e.g. 10-dimensional, representations, the corresponding SUSY-breaking masses need not be the 
same, i.e. the Higgs masses need not be universal [III [53, [28, 29]. In the scalar part of the general MSSM Lagrangian, 

£ so ft C - (ua u QH U -da<i QH d - ea e LH d + h.c.J 

— rriQ Q — V* vc^ L — uvo^'u — dm^d — em§e^ 

- m 2 Hu H: t H u - m 2 Hd H* d H d - (bH u H d + h.c.) , (25) 

the trilinear scalar couplings ai and the SUSY-breaking scalar masses m? are still unified to Aq and too at the 
GUT scale, but the SUSY-breaking Higgs mass parameters to# u and mu d are in general different. Models with 
non-universal Higgs masses (NUHM) are therefore defined by the parameters mo, TO1/2, Aq, tan/?, sgn(/x), to^, and 
rriH d - Note that in our NUHM models the gaugino masses remain unified to to^ and electroweak symmetry-breaking 
(EWSB) is still achieved radiatively, albeit through modified renormalization group equations (RGEs). This leads to 
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TABLE I: High-scale parameters together with the corresponding neutralino relic density, the contribution from top and bottom 
quark-antiquark final states to the annihilation cross section obtained with micrOMEGAs 2.1, and the mass eigenvalues of the 
lightest neutralino and the lightest stop for our three selected NUHM scenarios. 





m [GeV] M 2 [GeV] A [GeV] tan/3 sgn(^) 


m Hu [GeV] m„ d [GeV] 


^cdm^ 2 ti bb 


m-o (GeV) m h (GeV) 


I 

II 
III 


500 500 10 + 
620 580 10 + 
500 500 -1200 10 + 


1500 1000 
1020 1020 
1250 2290 


0.118 21.0% 64.0% 
0.118 51.0% 
0.113 93.4% 


207.2 606.4 
223.7 923.8 
200.7 259.3 



a more constrained parameter space in the too — mi/a plane as compared to mSUGRA models and to the fact that 
the minimum conditions of the tree-level Higgs potential 

-26 

sin/3 = — 5 5 -— r , (26) 

m 2 Hu + m 2 Hd +2/J, 2 ' y ' 



2~ ~ tan 2 (3 - 1 

m 2 _ 2 , 2 , o,,2 



/i , and (27) 



mjj u + m 2 Hd + 2fi 2 (28) 

allow to replace the parameter b by tan (3, but not to determine the superpotential parameter \i and the pseudo- 
scalar Higgs boson mass vtla as functions of too and m 1 / 2 , as was the case in mSUGRA. It is possible to replace the 
free parameters to# u and TO# d by the low-scale parameters (i and m^- Although we don't make use of this fact in 
our analysis, it can be useful to consider /i and tjia as free parameters of the model, since they strongly influence 
the annihilation of neutralinos into quarks and especially the relative weight of each channel. Having as a free 
parameter means that one can find scenarios, where m 2 A ~ 2m 2 , and the Higgs-boson exchange dominates the cross 
section even for smaller values of tan (3. In mSUGRA, such a scenario was only allowed for values of tan/3 > 40. We 
study this case by choosing the benchmark point I given in Tab. HI Furthermore, the higgsino parameter /i influences 
the higgsino component of the neutralino. By making it larger, one can enhance the contributions from the ^°-boson 
exchange. This can be achieved, as discussed in Ref. (26j . by starting with large positive values of and mu d at 
the GUT scale. By virtue of the RGE evolution, m 2 H is driven to small negative values, while Tn 2 H remains positive, 
which results in a small value of /i as given by Eq. |27[) . This in turn gives rise to a larger higgsino fraction of the 
neutralino and enhances the coupling to Z°-bosons. Such a scenario corresponds to our point II in Tab. |TJ Finally, 
one can make use of the NUHM RGEs to reduce the values of m^- and m5L This can be obtained by choosing a large 

and negative difference m 2 H ^ — m 2 Hd . On top of that, by choosing A to be large and negative, one induces a larger 
splitting in the third-generation sfermion masses, leading to a scenario (point III in Tab. [TJ) with a small top squark 
mass and enhanced t- and it-channel exchanges of top squarks. 



B. SUSY models without gaugino mass unification 

In the gaugino sector, the situation is somewhat similar to the one in Sec. IIII Al The breaking of the SO(10) 
symmetry can proceed via a step involving its SU(5) subgroup. The SU(5) later breaks into the Standard Model 
gauge groups SU(3)xSU[2)xU(l) [3(| and also determines the properties of the SUSY-breaking mechanism. As 
pointed out in Refs. [U, [22J, the breaking of supersymmetry itself is induced by an i^-term, while the gaugino masses 
are generated through a chiral superfield $, whose auxiliary component F$ acquires a vacuum expectation value. 
Assuming that the gauginos belong to the adjoint representation 24 of SU(5), the fields $ and F$ can in principle 
belong to any of the irreducible representations appearing in the symmetric product 

(24<g> 24) sym = 1© 24 ©75 ©200, (29) 

or any linear combination thereof. The relations between the gaugino masses Mi at the unification scale are given by 
the embedding coefficients of the Standard Model groups in SU(5). Note that these possibilities are all compatible 
with gauge coupling unification, but only the case where the SUSY-breaking field F$ is taken to be a pure singlet (1) 
leads to the gaugino mass universality featured by the mSUGRA model. 

The soft SUSY-breaking Lagrangian contains mass terms for the £?-ino, W-ino, and gluino, 

£ 80 ft c - ^M 1 BB + M 2 WW + M 3 gg + h.c^. (30) 
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TABLE II: High-scale parameters together with the corresponding neutralino relic density, the contribution from top quark- 
antiquark final states to the annihilation cross section obtained with micrOMEGAs 2 . 1, and the mass eigenvalues of the lightest 



neutralino and the lightest stop for our two selected non- 


universal gaugino mass 


scenarios. 






m [GeV] M 2 [GeV] A [GeV] tan/3 s 




X\ x 3 




tt 


m-o (GeV) 


m h (GeV) 


IV 


320 700 -350 10 


+ 


2/3 1/3 


0.114 


79.2% 


183.4 


281.9 


V 


1500 600 10 


+ 


1 4/9 


0.104 


50.4% 


235.6 


939.0 



As argued above, the values of Mi, M2, and M3 at the unification scale can be considered as independent parameters. 
Here, we adopt a commonly used parametrization and introduce the dimensionless parameters 

Mi M 3 
Xl = — and x 3 = — , (31) 

which will be used together with the lU-ino mass parameter M2 to describe the gaugino sector. The case x\ = x% = 1 
reproduces the mSUGRA model with the five parameters too, mi 12, A$, tan/3, and sgn(/z) at the unification scale. 

Models with non-universal gaugino masses have been shown to favor annihilationprocesses where neutralino annihi- 
lation into quarks is mediated by a Z°-boson exchange and squark exchanges [3l|, H2, HH . It has also been shown that 
the key parameter in this context is the gluino mass M3, since it influences practically all sectors of the low-energy 
mass spectrum through the renormalization group evolution. A decrease in M3 induces a decrease of the mass squared 
m 2 H , which through electroweak symmetry breaking conditions induces a decrease in the higgsino mass parameter 
fi as well. This in turn increases the higgsino fraction of the neutralino and lowers the pseudo-scalar Higgs mass 
mji. Moreover, having M 3 independent of the other gaugino mass parameters, it is straightforward to obtain lighter 
squarks, in particular the scalar tops. This effect can still be enhanced by decreasing the scalar mass parameter Too 
and by adjusting the trilinear coupling Ao, which influences the squark mass splitting. With scalar tops being light 
(even becoming eventually the next-to-lightest SUSY particle) , the squark exchange dominates the cross section for a 
B-ino-like neutralino, where a low value of tan/3 suppresses the Higgs-boson exchange (our point IV in Tab. [H}- As 
the higgsino fraction of the lightest neutralino increases, the squark exchange is enhanced due to the large Yukawa 
couplings. In the case of a large higgsino fraction, also the Z°-boson exchange becomes important and can take over 
if the lightest stop is not too close in mass to the lightest neutralino. We analyze the consequences of such a scenario 
by choosing the point V in Tab. HU 



IV. NUMERICAL RESULTS AND DISCUSSION 

Starting from the high-scale parameters, we use the public computer program SPheno 2.2.3 [3~i| for the numerical 
evaluation of the renormalization group running in order to obtain the SUSY-breaking parameters at the electroweak 
scale. The relic density of the neutralino is then evaluated numerically using the public code micrOMEGAs 2. 1 [Io| . 
where we have included our calculation of the neutralino annihilation cross section into third-generation quarks as 
discussed in Sec. [TTJ For the Standard Model input parameters, such as masses and couplings, we refer the reader to 
Ref. [35}, except for the value of the top quark pole mass, m t0 p — 172.4 GeV, which has been taken from Ref. (l2l |. 

We have chosen five typical parameter points shown in Tabs. HI and HI1 that have a dominant neutralino annihilation 
into top or bottom quarks through Higgs-boson, iT -boson or squark exchanges and whose relic density lies reasonably 
close to the WMAP range of Eq. JTJ). In all scenarios we only consider relatively low values of tan/3 = 10. In Sec. II V Al 
we analyze the first three parameter points (Tab. [TJ where we make use of the possible non-universality of the Higgs- 
boson masses to construct scenarios with cross sections dominated by Higgs-boson, Z°-boson, and squark exchanges, 
respectively. In Sec. IIVB1 we investigate the point IV, which corresponds to one of the scenarios without gaugino mass 
unification discussed in Ref. [3l[, and the point V, which is motivated by one of the "compressed SUSY" scenarios 
proposed in Ref. [32| . Our chosen parameter points also satisfy electroweak precision and low-energy constraints such 
as the measurements of the p-parameter, the anomalous magnetic moment of the muon, and the branching ratio 
of the decay b — > sj. Due to the large experimental error on the measurement of Ap [35j, however, only regions 
featuring very high SUSY masses are excluded at the 2cr-level, so that this constraint does not affect our analysis. We 
have taken into account the anomalous magnetic moment of the muon by taking the higgsino mass parameter p, > 0. 
Negative values are disfavored, since they increase the gap between the recent experimental value and the Standard 
Model prediction [36| . The most stringent constraint is given by the inclusive branching ratio of the decay b — > S7. 
Recent experimental measurements from BaBar, Belle, and CLEO lead to the combined value (3?| 

BR(6^s 7 ) = (3.52 ± 0.25) • KT 4 . (32) 
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The theoretical prediction of the SUSY contribution is particularly sensitive to the masses of the chargino, the charged 
Higgs boson, and the lightest scalar quark. We have verified that our parameter points lead to values that lie within 
2a of the above limit using the public codes FeynHiggs 2.6.5 [HI and SusyBSG 1.3 [33]. Note that also the direct 
mass limits from collider searches are fulfilled [35j . 

Although the focus of the present paper is the improvement of the annihilation cross section through SUSY-QCD 
corrections, one should keep in mind that there are also other sources of uncertainties in the calculation of the relic 
density. From the particle physics side, it is well known that differences in the low-energy mass spectrum may occur 
when using different spectrum generators. This, in consequence, may induce a sizable difference in the prediction 
of the dark matter relic density or other observables Q. Note that, in this context, also the numerical value of the 
Standard Model parameters, in particular the top-quark mass, can influence the favored regions of parameter space, 
in particular in the case of dominant annihilation into top-quark final states. Corrections to the annihilation cross 
section are also induced by the electroweak interaction Q , but these are generally smaller than the strong corrections 
and beyond the scope of this work. Concerning cosmology, it has been shown that modifications of the standard 
cosmological model may also affect the prediction of the dark matter relic density. Including, e.g., an energy content 
such as quintessence or an effective dark density due to extra dimensions modifies the expansion rate or the entropy 
density of the early universe (see e.g. and thus enters into the calculation of the relic density. 

In the left panels of Figs. OH [Jj and [lOl we show the cross section for the annihilation of a neutralino pair into 
a bottom and/or top quark-antiquark pair as a function of the center-of- momentum energy p CTa calculated at the 
tree-level with DR Yukawa couplings (solid lines) as well as the leading contributions from the individual annihilation 
channels (dashed lines) and their interferences (dotted lines). Note that some of the latter have been multiplied by 
(— 1) in order to fit into the logarithmic plot. We also show, in arbitrary units, the thermal velocity distribution 
function, involved in the calculation of the thermal average evaluated at the freeze-out temperature (shaded areas). 

In the right panels of Figs. l5l 171 and [TUl we show again the total annihilation cross section into bottom and/or top 
quark-antiquark pairs at different levels of precision, i.e. at leading order with DR couplings (dash-dotted lines), in the 
approximation included in the micrOMEGAs code (dashed lines), which uses effective couplings to absorb the leading 
loop effects, and with our full one- loop QCD and SUSY-QCD corrections (solid lines). The shaded area indicates 
again the thermal velocity distribution of the neutralinos at the freeze-out temperature in arbitrary units. 

In order to generalize our results, the remaining Figs. [6l [8J El HTl and [T2l show scans in various two parameter planes 
around our parameter points in Tabs. HI and HT1 and display the contours allowed by WMAP calculated at tree- level 
(green), with the approximation implemented in micrOMEGAs (red), and with our full SUSY QCD corrections (blue). 
We also show dependence of the relic density on various physical parameters, e.g. the mass of the lightest neutralino. 

A. Relic density in models with non-universal Higgs masses 

We begin our detailed numerical discussion by analyzing a scenario, where the annihilation cross section into quarks 
is dominated by an exchange of a Higgs boson. This scenario corresponds to our parameter point I from Tab. HI Here 
the heavy CP-even and the CP-odd Higgs boson resonances coincide with both Higgs boson masses at 450.3 GeV. 
Due to the dominance of the Higgs-boson exchange, the neutralino annihilates predominantly into bottom quarks. 
This is a consequence of the fact that, although the bottom quarks have a much smaller mass compared to the top 
quarks, their couplings to the Higgs bosons are enhanced by tan/3 or cos a. As opposed to mSUGRA, the neutralino 
in this scenario is rather heavy, which allows for an effective annihilation into top quarks as well. In the bottom 
panels of Fig. [SJ a small discontinuity indicates the place, where the top quark final state starts to contribute. These 
circumstances result in a mixture of top and bottom quark final states, which cannot be reached in mSUGRA for 
such a low value of tan/3. A few general features can be observed in Figs. [5] and [5] First, on-resonance annihilation 
of neutralinos would reduce the relic density too much, so that regions allowed by Eq. (fT]) sit on each side of the 
resonance peak. In Fig. O one can clearly see that the resonance is not aligned with the kinematic region where the 
biggest contribution to the relic density comes from. The top panels in Fig. [6] display the bands where quark final 
states dominate the annihilation. In addition, we show a line that corresponds to the position of the Higgs resonance 
peak. The WMAP allowed regions are situated on each side of the resonance where the width of the these regions 
reflects the steepness of the peak. The SUSY-QCD corrections to the relic density calculated here are a combination 
of corrections to the processes with either bottom or top quarks final states. The corrections to the annihilation into 
bottom quarks are essentially confined to the Higgs-quark-antiquark vertex. The bulk of the correction comes from 
the gluon exchange and from the SUSY corrections which become large at large tan (3 and have to be resummed. As 
we take a moderate tan/3 = 10 even the non-resummable corrections play a role here. These corrections amount to the 
difference between our full SUSY-QCD result and the effective coupling approximation implemented in micrOMEGAs 
in this scenario. The situation with the top quarks is more complicated. As can be seen from Fig. [5j the squark 
exchanges in the t- and it-channels and their interference with the Higgs-boson exchange become now also sizable, 
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FIG. 5: The contributions of the different diagrams to the annihilation cross section of a neutralino pair into bottom (above) and 
top (below) quark-antiquark pairs (left) and the effect of the radiative corrections on the annihilation cross sections (right) as a 
function of the center-of-momentum energy p cm for our parameter point I. The shaded area indicates the velocity distribution 
of the neutralino at the freeze-out temperature in arbitrary units. 



which renders also the corrections to squark-quark-neutralino vertex and to the squark propagator important. Overall 
for this scenario, the effect of the full SUSY-QCD corrections amounts to a 20% shift in the prediction of the relic 
density with respect to the calculation implemented in micrOMEGAs, as can be seen in the bottom panels of Fig. [6] In 
consequence, the cosmologically favored regions of parameter space are shifted away from the position of the Higgs 
pole in order to compensate the larger annihilation cross-section. These shifts amount to up to 5 GeV for the common 
scalar and gaugino masses and a few GeV for the masses of the Higgs doublets at the unification scale (see Fig. [|5]). 
This is about two times larger than the experimental precision, and therefore distinct bands appear. 

The remaining analyzed scenarios were chosen so as to investigate effects not related to Higgs-boson exchanges. 
The Higgs-boson masses are very heavy in these scenarios (m H 0,771^ > 1100 GeV) and suppress the importance 
of the Higgs-boson s-channel exchange. The parameter point II in Tab. Q] was chosen as described in Sec. IIII Al so 
that the Z°-boson exchange is enhanced. The only viable way to do this is to increase the higgsino fraction of the 
neutralino, which in turn increases the coupling of neutralinos to the Z° boson. Another possibility would be to sit 
on the Z°-boson resonance, which would, however, lead to very small neutralino and chargino masses that are already 
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FIG. 6: Top: Cosmologically favored regions in the mo-m^-plane (top left) and rn,H u -Tn.H d -plane (top right) for scans around 
our parameter point I. We show the regions that satisfy the constraint from Eq. |T} for our tree-level calculation (green), the 
calculation implemented in micrOMEGAs (red), and our calculation including the full SUSY-QCD corrections (blue). We also 
indicate the contributions from quark-antiquark final states to the total annihilation cross section by isolines. Bottom: The 
prediction of the neutralino relic density QcdmIi 2 including the tree-level (green dash-dotted) cross section, the approximation 
included in micrOMEGAs (red dashed), and the full one-loop SUSY-QCD corrected cross-section (blue solid) as a function of the 
gaugino mass parameter mi/2 (left) and as a function of rn,H u parameter (right). The shaded area indicates the favored region 
of Eq. gj. 



ruled out by the LEP experiments. The parameter point III in Tab. [T] was picked to use the Higgs potential high-scale 
parameters to drive down the squark masses via the renormalization group evolution. On top of that, by choosing 
Aq = —1200 GeV we induced a large mixing of the third generation squarks. The mixing is biggest for scalar tops 
increasing the contribution of their exchange. A feature common to both scenarios is a distinct destructive interference 
between the squark and the Z° boson which links these two scenarios (see Fig. [7]). This fact makes corrections to 
both Z°-boson and squark exchange important in each of the two scenarios. In the case of the parameter point II, the 
corrections are about 20% in the cross section as compared to the cross section in micrOMEGAs. This is not reflected 
in all regions of the WMAP allowed regions in Fig. [5J since the contribution of the quark-antiquark final state falls 
to only 60%-40%, as we lower the higgsino parameter [i. This decreases the masses of the lightest chargino and of 
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FIG. 7: The contributions of the different diagrams to the annihilation cross section of a neutralino pair into top quark- 
antiquark pairs (left) and the effect of the radiative corrections on the annihilation cross section (right) as a function of the 
center-of-momentum energy p cm for our parameter points II (above) and III (below). The shaded area indicates the velocity 
distribution of the neutralino at the freeze-out temperature in arbitrary units. 



the second-lightest neutralino, which, through their t-channel exchange, increase the annihilation cross sections into 
W + W~ and final states. Nevertheless, the full SUSY QCD corrections shift the contour in the rriA — £t plane 

by 5 GeV in ii and in all instances shift the contour by more than the current experimental precision. The effect of 
the corrections to the cross section is not screened by other final states in the case of the scenario with a dominant 
squark exchange (point III in Tab.H|. The top quark final state accounts in certain regions for more than 90% of the 
annihilation cross section, and the mass of the scalar top is not light enough to allow for efficient co-annihilations. The 
hyperbolic shape of the WMAP allowed regions in the — m H d plane is governed by the Higgs-mass parameter 
combination m 2 H — rn^ d . The effect of the corrections on the cross section is about 20% and causes a shift of the 
preferred value of the pseudo-scalar Higgs boson mass by about 50 GeV. 
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FIG. 8: Top and Bottom left: Scans around our parameter point II in the mo-mx/2-plane, mfr^-m^-plane, and m^-/i-plane. 
Bottom right: The prediction of the neutralino relic density ficDM^ 2 as a function of the high-scale Higgs-mass parameter 
itih u ■ The labels are the same as in Fig. [6] 

B. Relic density in models without gaugino mass unification 

The relevant masses of our points III and IV featuring very light stops are very similar, as can be seen in Tabs. U 
andflj] The same holds for the annihilation cross sections. For the benchmark point IV, the exchange of a stop in the 
t- or it-channel is therefore favored, which is well visible in the top left panel of Fig. [10] The contribution from the 
Z°-boson exchange is here one order of magnitude lower than the one from squark exchange. Again, an important 
role is played by the squark Z°-boson interference effect, which is here even stronger than in the case of point III. The 
differences are the Higgs boson masses, which were more than 2000 GeV for point III, but are m H o ~ 622 GeV 
here (see Fig. [TO]). The mass difference can be traced back to the mechanism by which we lowered the squark masses. 
For point III, we had to induce a big difference between the high-scale Higgs potential parameters, which drove the 
Higgs-boson masses higher, whereas for point IV we changed freely the gluino mass parameter influencing the RGEs 
for the squark masses, which has not such a big effect on the masses of the Higgs bosons. Due to the large contribution 
from the exchange of the light squark, the annihilation into top quarks accounts for up to 80% of the total annihilation 
cross section (see Fig.Ql]). The subdominant final states here are again W + W~ and Z°Z°, but also pairs of bottom 
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FIG. 9: Top and Bottom left: Scans around our parameter point III in the mo-mi/2-plane, m^-ma^-plane, and mA-/i-plane. 
Bottom right: The prediction of the neutralino relic density Ocdm^ 2 as a function of the high-scale Higgs-mass parameter 
m,H u ■ The labels are the same as in Figs. [6] and [8] 



quarks or leptons. Note that also here the mass difference between the lightest neutralino and the lightest stop is not 
small enough for efficient co-annihilations. One small, but important difference is that the destructive interference is 
stronger for point IV than for point III. This leads to negative correction effects when including the full SUSY QCD 
corrections (when compared to the cross section of micrOMEGAs). The reason is that for this particular point the 
corrections to the squark exchange are effectively taken into account by including DR Yukawa couplings and the only 
difference between our full calculation and the micrOMEGAs approximation are the corrections to the interference which 
push the full result down (as seen in the top right part of Fig. 1101) . Let us now turn to the impact of the corrections 
in this "compressed SUSY" scenario. The corresponding favored regions in the mo-Mi and the tuq-xj, planes, the 
projection on the m^-m^ plane as well as the prediction of the neutralino relic density as a function of the gaugino 
mass parameter Mi around our parameter point IV are shown in Fig. 1111 Again, in the regions where annihilation 
into a top quark-antiquark pair is kinematically allowed, the effect of the SUSY-QCD corrections is sizable, resulting 
in an important shift of the favored regions in the parameter space. The preferred region of parameter space is shifted 
to higher values of the scalar mass parameter m and consequently to higher stop masses . As already discussed 
above, a heavier stop mass compensates the increasing effect of the additional loop diagrams on the annihilation 
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FIG. 10: The contributions of the different diagrams to the annihilation cross section of a neutralino pair into top quark- 
antiquark pairs (left) and the effect of the radiative corrections on the annihilation cross section (right) as a function of the 
center-of-momentum energy p cul for our parameter points IV and V. The shaded area indicates the velocity distribution of the 
neutralino at the freeze-out temperature in arbitrary units. 



cross-section. In the right bottom panel of Fig. [TT] one sees again that the prediction of the relic density is decreased 
by the order of 15% with respect to the tree-level calculation and by about 10% with respect to the approximation 
implemented in micrOMEGAs. 

For our point V with a dominant Z°-boson exchange, the correction accounts for about 20% of the micrOMEGAs 
cross-section. Although micrOMGEAs does not include any correction for the exchange of a Z°-boson, the corresponding 
curve is approximately 20% over the tree-level prediction, which can be explained by the presence of effective couplings 
for the sub-leading squark-exchange. The difference between the approximation included in micrOMEGAs and our full 
one-loop calculation originates from the supplementary corrections, especially for the exchange of the Z°-boson, shown 
in Figs. [Hand [3]and discussed in Sec. [TTJ We study again the influence of our corrections on the regions of parameter 
space that are favored with respect to the WMAP limits of Eq. (TTJ). The upper left panel of Fig. [T2l shows these 
regions in the m^-M-i plane for our scenario V and fixed values of x\ and x% as given in Tab. [TTJ In the lower left 
panel, we show the same contours projected on the corresponding plane of the physical neutralino and stop masses 
rrix and m^. In each plot, we also indicate the isocontours corresponding to a contribution from top quark- antiquark 
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FIG. 11: Top and bottom left: Cosmologically favored regions in the m,o-M2-plane (top left), mj-mj-plane (bottom left), and 
rno-a;3-plane (top right) for scans around our parameter point IV. We show the regions that satisfy the constraint from Eq. ([1]) 
for our tree-level calculation (green), the calculation implemented in micrOMEGAs (red), and our calculation including the full 
SUSY-QCD corrections (blue). We also indicate the contributions from top-quark final states to the total annihilation cross 
section by isolines. Excluded regions due to unphysical solutions of the renormalization group equations are shown in gray. 
Bottom right: The prediction of the neutralino relic density QcoMh 2 including the tree-level (green dash-dotted) cross section, 
the approximation included in micrOMEGAs (red dashed), and the full one-loop SUSY-QCD corrected cross-section (blue solid) 
as a function of the high-scale gaugino mass parameter M2 for fixed mo = 300 GeV. The shaded area indicates the favored 
region of Eq. ([1]). 



final states of 50%, 30%, and 10% to the total neutralino annihilation cross-section. The points lying in the grey 
shaded areas do not allow for physical solutions of the renormalization group equations. The correction to the relic 
density is reduced from the 20% it was at the cross section level (as compared with the micrOMEGAs cross section), 
because the quark-antiquark final state constitutes only about 50% of the annihilation cross section. The remaining 
contributions include mostly W + W~ and Z°Z° final states. Nevertheless, as it was also the case for the mSUGRA 
scenarios analyzed in Ref. [9|], the impact of the one- loop SUSY-QCD corrections is larger than the experimental 
uncertainty. Therefore distinct bands are observed in wide regions of both the m^-M-i and in the m^-m^ planes. 
Note that for M2 < 450 GeV, the annihilation of a neutralino pair into top quarks is kinematically forbidden. The 
dominating channels for that region are mainly the annihilation into combinations of the gauge bosons and the light 
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FIG. 12: Top and Bottom left: Scans around our parameter point V in the mo-M2-plane (top left), m^-m^-plane (bottom left), 
and mo-a;3-plane (top right). Bottom right: The prediction of the neutralino relic density QcDwih 2 as a function of the stop 
mass mj for a fixed neutralino mass = 235.6 GeV. The labels are the same as in Fig. 1111 

Higgs boson. In the m^-x^ plane, shown in the right top panel of Fig. [T2l the WMAP-favored points are confined to 
rather narrow bands corresponding to the different levels of included corrections. This underlines that M% is one of the 
key parameters to which phenomenology is very sensitive. It is interesting that for lower values of x% their positions 
are almost independent of mo, while for higher X3 the dependence becomes stronger. The shift from the tree-level 
prediction to our full one- loop result is in the direction of higher gluino masses M3 = 23 A/2- This is explained by 
the fact that increasing the gluino mass implies an increase in the squark masses, which in turn implies a decrease of 
the annihilation cross-section. This effect is then combined with the increase of the cross-section due to the one-loop 
corrections discussed above, so that the relic density remains in the narrow range which is favored by cosmological 
data. Finally, in the right bottom panel of Fig. [12l we show the prediction for the neutralino relic density Sl^h as 
a function of the lightest stop mass . The neutralino mass has been fixed to the value = 235.6 GeV of our 
point V. The graph corresponds to a cut through the m^-m^ plane shown in the lower left panel of Fig. 1121 The 
favored region of Eq. (TTJ) is indicated by a shaded area. Since the SUSY-QCD corrections increase the annihilation 
cross section by about 50%, the prediction for the neutralino relic density is reduced by about the same amount. 
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V. CONCLUSIONS 

The theoretical calculation of the dark matter relic density is an interesting tool to obtain rather stringent constraints 
on the MSSM parameter space, both at the electroweak and at the grand unification scale. Cosmological precision 
measurements therefore play an important role in the extraction of SUSY mass parameters from experimental data. 
In times of increasing experimental sensitivity, in particular for the cosmological parameters such as the relic density 
of cold dark matter, it is essential to increase the accuracy of the theoretical calculation. Consequently, higher- 
order corrections to the dark matter annihilation cross section become important, since they enter directly into the 
calculation of the relic density. 

We have presented here the full analytic details of our calculation of the 0(a s ) corrections to the annihilation 
of a neutralino pair into a massive quark- antiquark pair, i.e. of the one-loop and real gluon emission corrections. 
Annihilation processes into heavy quarks have been shown to be important in large cosmologically allowed regions of 
SUSY GUT models without scalar or gaugino mass unification, and the effects of our SUSY-QCD corrections have 
been investigated numerically for these two classes of models. We have identified regions of parameter space, which 
correspond to current cosmological limits on the dark matter relic density and which feature important annihilations 
of the neutralino pair into a top quark-antiquark pair through the exchange of a Z°-boson or a scalar quark. Higgs 
boson resonances were shown to be suppressed in those regions. For five selected parameter points, the effects of the 
corrections were shown to be sizable, since they enhanced the annihilation cross-section by typically 20 and up to 
50% with respect to the tree-level calculation. 

As a consequence, the theoretical prediction of the neutralino relic density is also affected by the contributions at 
the one-loop level. Since the cross section is increased by typically 20 and up to 50%, the relic density is reduced by 
about the same amount. We have shown that the impact of our corrections is more important than the uncertainty on 
the observational limits at the 2a confidence level. It is therefore essential to take these corrections into account when 
analyzing cosmological data with the goal of extracting SUSY mass parameters, which may be shifted by typically 5 
GeV and up to 50 GeV, and of determining the favored regions of the SUSY parameter space. 
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APPENDIX A: NOTATION AND COUPLINGS 

We follow closely the notation of the couplings defined in Ref. and the conventions used in Ref. [lj, [lj| and 
begin by listing all necessary couplings of the Z°-boson. The couplings of fermions to the ^°-boson in the MSSM are 
identical to those in the SM, and the Lagrangian is 

C = -^zlfr(C[P L + C f R P R )f, (Al) 

Cw 

where g is the weak coupling constant, = Ij L — CfS^, C^ R = — e/s^,, s w (c w ) is the (co-)sine of the weak mixing 
angle 6 W , and Ij L and e/ are the weak isospin and electric charge of the fermion /. The Lagrangian for the interaction 
of the Z°-boson with two neutralinos is given by 

c = 2^ Z ° K LPi + °'^ Pr ) % ' (A2) 

where 

= — ^ZisZjs + T^ZaZji = — O io R (A3) 
depend bilinear ly on the neutralino mixing matrix Z . The Lagrangian of the Z°-boson coupling to two sfermions 
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is proportional to 



J — nf v>f r?f 



which depends bilinearly on the sfermion mixing matrix 
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We continue with the couplings containing neutral Higgs bosons, where we use the notation H% — {h°, H°, A , G }. 
t/t stands for an up-type (s)fermion and b/b for a down-type one. For the neutral Higgs-boson-fermion-fermion 
couplings, the interaction Lagrangian reads 
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(A7) 



(A8) 



Here, h t and h b are the Yukawa couplings 



h t = 



\/2»Ti lv sin/3 



hh = 



V2m w cos {3 



The interaction Lagrangian for neutral Higgs bosons and neutralinos is given by 

2 4 
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where dk and take the values 

dk = {— cos a, — sina, cos/3, sin/3} , = {— sina, cos a, — sin /3, cos/3} 

Following Ref. [40] , the neutral Higgs-boson-sfermion-sfermion couplings can be written as 
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where the left-right couplings G LR k for third-generation up- and down-type sfermions are 
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^=(A t c a + ns a ) 



'\/2h t m t c a + g z m z e t sl,,s a+ p 
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Gb 
LR.l 



^LR.2 



^LR,3 



V2 

' u (A b s a + \xc a ) V2h b m b s a + gzm z ebsl,-s a+ p 



V2h b m b s a + g z m z (ll L -e b s 2 w )s a+ p 
111. 

with a — > a — tt/2 , 




%(A b s a + ^c a ) 



nf 

u LR,l 



G 



LR,3 



-V2h t 
-V2h. 



%(A t C/3 +jUS /3 ) 









G 



f 

LR,4 



G 



f 

LR,3 



(A b sp + \ic$) 
with J3 -> P - tt/2 



|(ylb s/j + pep) 




Here, we have used the abbreviations s x = sin x and = cos x and a denotes the mixing angle of the {h , iJ°}-system. 
For the neutralino-sfcrmion-fcrmion couplings, the Lagrangian reads 



C = f ("IPr + b{ k P L ) xl fi + Xk (af k PL + b{ k P R ) If* 

with the coupling matrices 

4 = I'f ■ //J/,":, • b{ k = hfZ^+gf^R^ 

and 

f[ k =V2((e f -I 3 f L )t<m6 w Z kl +I 3 f L Z k2 ) , f f Rk = -V2e f t<m8 w Z kl . 
The coupling of gluinos to sfermions and fermions can be derived from the Lagrangian 

C = -V2g s T« \f s {R lL P R - R lR P L ) g a f i<t + ~g a (R lL P L - R lR P R ) fj* 



(A12) 

(A13) 
(A14) 

(A15) 



where g s = y/-iira s is the strong coupling constant and R is the sfermion mixing matrix defined above. Finally, the 
Lagrangian corresponding to the couplings of the gluon to fermions and sfermions is 



C = g s T»g a J s ^f t igJ^g-f-jy'fjj. 
where the T a represent the usual SU(3) color matrices. 



(A16) 



APPENDIX B: VIRTUAL CORRECTIONS 



Here, we give explicit forms of the next-to-leading order amplitudes A4i_i 00 p mentioned in SecHU We use generic 
amplitudes with generic couplings, which wc then specify using the couplings given in App.fX] Note that all kinds of 
indices are to be summed over even if not explicitly stated. We leave out the results for self-energies, which can be 
found in 19| using the same conventions and couplings. 



Vertex corrections 



The contribution of the vertex corrections depicted in Fig. [2] to the matrix element A^i_i 00 p can be written in 
terms of six loop diagrams. Denoting the momenta of the incoming neutralinos p\ and p 2 and those of the outgoing 
quark and anti-quark as k\ and k 2 , the next-to-leading corrections to the s-channel exchange of the J?°-boson can be 
parameterized as 



■Mi-ioop = V(p2h"{A L P L + A R P R )u( Pl ) 



(Bl) 



t(fc 1 ){ 7 I '(BiP L + B R P R ) + k\(B 2 L P L + B 2 R P R ) + kl{B\P L + B 3 R P R )}v(k 2 ) 



22 



q+ki 



V„ 



q-k 2 



V 3 *' 



V * k 2 



v* fa 



V 3 :7 P (93 L ^+33^fl) 
Vi :7 P (54 L PL+5fPfl) 




" 1 



Vo 



ff0 (2g + fci -k 2 ) v 
(93 Pl + 9iP R ) 
{gkP L + g?P R ) 



v 4 



FIG. 13: Vertex corrections to the Z°-boson coupling to quarks, where the arrows indicate the flow of the four-momenta and 
q is the independent loop momentum. 



where Al , Ar are the chiral couplings of two neutralinos to the Z°-boson 



A t 



-A f 



2cw 



and B l L , B l R are general form factors, which will now be given for the loops with a gluon and a gluino exchanj 
form factors for the loop diagram containing a gluon (on the left side in Fig. 113(1 are 



n l 1 L L 

B L = ^g 39i 



B 2 T = 



Bi = 



Air 2 



, 9o L (B - 2C 00 - s(C + Ci + C 2 ) + m) (2C + 3Ci + 3C 2 )) + g$m) {d + C 2 ) 
g^g^m f f^di + g£(C 12 + C 2 ) 



,L„L, 



4tt 2 



Br = Biting*), i = 1,2,3. 

Using the notation from App. the generic couplings in this case are 



go — °li 

Cw 



R_ J_ n f 

cw 



L R L R 
53 = 53 = 54 = 54 



... St 



The scalar loop integrals follow the definition given in Ref. [4l[ and their arguments are 

Ci(ki,-k 2 ; 0, m/,m/) , Bo( - fci - fc 2 ; m/, m/) , 

where m/ is the mass of the final-state fermion. 

The same form factors for the diagram with a gluino exchange (see Fig. [T3] right) are given as 



50 g^g a Coo , 



B\ = 



B\ 



1 

16tt 2 
1 



50 



■ 5o 



16tt 2 

Br = BKg^g^g^gK). 



5 4 L ( - 5 3 WO) + 2d) + <?f m/Cd + 2C U )) - &9%m f {(h + 2Ci 2 ) 

5 4 L ( - 53 i ^5(Co + 2C 2 ) + ^m/(Ci + 2C 12 )) + 5?53™/(d- + 2C 22 ) 
1 = 1,2,3, 
where the couplings are 

50 = -— 4' 93= V2g s T-R. l2l 

c w 

g£ = ~V2g s T: t R 3l , 

Here, the arguments of the loop integrals are 

Ci(k-L,—k 2 ;mg,m j. ,mj.) , B ( -ki -k 2 ;m^,m^) 



g* = -V2g s T? t R tl 



9? 



V2g s T°R j2 
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(B14) 
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FIG. 14: Vertex corrections to the ffj:-boson coupling to fermions, where the arrows indicate the flow of the four-momenta and 
q is the independent loop momentum. 



The generic structure of the amplitude, which parameterizes the corrections to the s-channel Higgs-boson exchanges, 
is much simpler and can be written as 



Mi — loop — 



k=l 



'tl, 



l{kx){Dl <k P L + D R<k P R )v(k 2 ) 



where Cl± and C Rik are the couplings of the two neutralinos to the Higgs boson 



Clm = Ci?.,fc — 



9 p o 

C L ,k = -C Rt k = | Ff lk 



for k = 1,2. 



for k = 3, 4 . 



(B15) 

(B16) 
(B17) 



There are only two form factors T)\ and D Rl which receive contributions from two diagrams (see Fig. 1 14[) . The 
contribution in the case of the gluon exchange is 



n l 1 LL 



gk(2B - m}C + (s + 3m})(C + d + C 2 )) - g^m}(Ci + C 2 ) 



where the couplings are 

9o=9^ = 4 {k = 1,2), 



9 L n 



-9? 



-4 fc=(3,4), g L = gK = gl=gK = -g s T* t . 



The loop integrals are identical to those defined in Eq. (|B8 
The same form factors for the gluino contribution are 



Ion' 1 



g^g^m- g Co + g^g^m s C x + g%gfm f C 2 



L„R„ 



Dr = Dl(g^g^g^gi), 



where the couplings are 



.90 = G{ 



ijk ' 



.93 

9l 



= V2g s T«R l2 



.93 



= -V2g.T^Ra , 



(B18) 
(B19) 

(B20) 

(B21) 
(B22) 

(B23) 



-s/2g s T: t R jx , gf = V2g 8 T^ t R j2 , 



and the loop integrals are defined in the same manner as in Eq. (iBTill . 

The remaining vertex corrections are those connected with the t- and w-channel squark exchanges, where we have 
defined t = (pi — k 2 ) 2 and u = (pi — fci) 2 . We give here just the results for the i-channel, since the it-channel can 
obtained by using the crossing symmetry. The generic amplitude is 



Mi- loop = Y,^ k i){ E l,iPL + E l R>i P R )u(p 2 ) 

i=l 
2 

+ J2 [u(ki)(E L ,iP L + E R ,iP R )u(p 2 ) 



t — mi 



v(pi)(F L ,iP L + F R<i P R )v(k 2 ) 
VfaKFLPL + FkiPRWha) 



(B24) 
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FIG. 15: Vertex corrections to the neutralino coupling to a fermion and a sfermion where the arrows indicate the flow of the 
four-momenta and q is the independent loop momentum. 



where -Ex/K.i and Fl/rj are tree-level couplings given as 

El.i = b{\ , — d-d > Fh,i = Eft.i 

The contribution of the gluon exchange diagram (see Fig. [15] left) is 
1 



R,i 



Er 



16tt 



1 [„L 



El 
E \ 

F L = E R> 
F r = E\ , 



2 9o9s 



Ei(9 



94 



g%(B + 2m)(C - C\) + (s + t)C x + m\(C 2 - C\)) - gfm^m f (C + C 2 + C\) 

R\ 



where the couplings are 

.90 = -QaTStSa , <? 3 L =gf = - 
and the scalar loop coefficients have the arguments 

Cii-h, -p 2 ;mf,0,mf.) 

The gluino contribution gives 
1 



9sT s a t 



9i 



9 4 = 4i 



Both - P2;0,m f ~.) . 



(B25) 

(B26) 

(B27) 
(B28) 
(B29) 

(B30) 
(B31) 



Eh 



Er 
Ft 
Fr 



9o93 94m S mfC + g£ g£ g^m- g m^(C + C 2 ) + g^g^mj(C + Ci) + g^g^m-gmfd (B32) 



L L „R 



L„R„L2 , 



167T 2 

„L R R„ 



+9o93 9?m imf {Co + C\ + C 2 ) + g£g^gX{B + mjC + mjCi + m\C 2 ) + <?o 9s 9% m s m ± C 2 



' • •'/ / :'.'/', • •'/':'.'/! 



E r , 
El, 



where the couplings are given by 



. 9o L = V2g s T: t R a . 



9 3 



54 = a f ji ■ 



9? 



9? 



-V^gsT^Ra , 
and the scalar loop functions are 

Ci(-ki, -p 2 ; mj. , m g ,m f ) 



V2g s T s a t R j2 , 
-V2g s T* t R n , gf = I 

Boih - p 2 ;m g ,m f ) . 



(B33) 
(B34) 
(B35) 



(B36) 



(B37) 



Box corrections 



At the one-loop level, the neutralino annihilation into quark-antiquark pairs receives box contributions arising from 
the exchange of a gluon or a gluino between the final state quarks (see Fig. [5]) . In order to express the corresponding 
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a%P L +b%P R 
bf k PL+al k P H 



V 2 



FIG. 16: Conventions and notations used for momenta, masses, and couplings in the calculation of the box diagram arising from 
the exchange of a gluon between the two final-state quarks in neutralino-pair annihilation. The arrows indicate the direction 
of defined four-momenta. The coupling strengths are defined in Sec. 1X1 



amplitudes in a rather generic way, we introduce the notations shown in Figs . [TBI and 1171 for the gluon and gluino box, 
respectively. The momenta of the incoming neutralinos are again labeled p\ and p 2 , those of the outgoing quarks are 
k\ and k 2 . The corresponding external masses are rrii (i — 1, ... ,4), while the masses of the internal particles are 
denoted Mj (i = 1, . . . ,4). The left- and right-handed coupling strengths, which correspond to the vertices appearing 
in the schematic diagram, are generically denoted by gf R (i = 1, . . . , 4). For the sake of compact expressions, we use 
in the following the abbreviations 

Vi = gtP L +gKp R , V, = gfP R + g*P L , (B38) 

where P^ and P R are the left- and right-handed chirality projectors, respectively. The independent loop-momentum 
q is defined to be the momentum of the particle having the mass Mi as shown in Figs. [TBI and IT71 The arising tensor 
integrals have been reduced to scalar integrals 

D {iiij} = D {hlj} (-k 1 ,-(p 1 +p 2 ),-p 2 ;M 1 ,M 2 ,M 3 ,M i ) (B39) 

for i, j = 0, ... ,3. For their definition and the applied convention we refer the reader to Ref. (4lj . 
In case of the gluon box, the amplitude can be written as 

3 2 

iM ( h 9 l = -u(k 1 )^u{p 2 )v{p 1 )^v(k 2 )D -X)5^(Wf - 1 u{p 2 )v(p 1 )^ k v(k 2 )D i 

i=i k=i 

3 

-u{ki)TlQu(p 2 )v{pi)Tl Q v{k2)Dw - ^2 u{ k i)Fij u (P2)v(px)FijV{k 2 )Dij, (B40) 

where the form factors are given by 

Tl = m l7 ^V 3 V 2 

T\ = M 3 ^V 3 V 2 - 2ti{V 3 V 2 
H = -m l7 ^V 3 V 2 
T\ = -7" ^V 3 V 2 - M 2 7 M V 3 V 2 
Tl = m l7 ^V 3 V 2 
Tl = -M 2 ^V 3 V 2 
Tl = m l7 fV 3 V 2 

Ho = VWV2 

T\ x = 2A£ V 3 V 2 - M 3l »V 3 V 2 
T\ 2 = 2fcf V 3 V 2 - M 37 ^V 3 V 2 
T\ 3 = 2k1 V 3 V 2 - M 3l »V 3 V 2 



Tl = m 3 Vi>V 4 + MxVx^Vi- ^ 2 Vi>V 4 

Tl = m 3 Vi 7p V 4 + M x Vi>V 4 - ^ 2 Vi 7m V 4 

T\ = #iVi 7ai V 4 

Tl = m 3 Vi 7p V 4 + Mi Vi 7jU V 4 - ^ 2 Vi 7m V 4 

Tl = M 1 V l7 , I V 4 - ^ 2 V i7m V 4 

Tl = m 3 Vi 7M V 4 + MiVi 7/1 V 4 - ^ 2 Vi 7m V 4 (B41) 

J| = ^ 2 Vi 7m V 4 

^00 = Vi> 7ai V 4 

T\ x = hVi"f^4 

T\ 2 = -MiVi 7ai V 4 + ?* 2 Vi 7/ .V 4 

^13 = &Vi 7 uV 4 



2G 




P, _L f,<2 



a%P L +b%P R 
V 3 : ~V2g s T^R jL P L - fljafii) 
14 : -V2g s T? t (RuiPl - R iL Pr) 



FIG. 17: Same as Fig. [16] for the exchange of a gluino between the two final-state quarks. 
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(B42) 



Using the same notation, the amplitude of the gluino box is given by 

3 2 



iM 



(§) 

box 



= 1 fc=l 



+ v(p 2 )^ 1 u(pi)u(fci)^"o ?;(fc 2 )Doo 



X! v (P2)^ijU{pi)u(ki)J 7 f j v(k 2 )D i 

i,j=l 



(B43) 



with the form factors 
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M 2 V 2 Vi 


n 


= m 2 V 3 V 4 + M3V3V4 
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= M3V3V4 + 


m 2 V 3 V 4 




= M2V2V1 + 
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